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Abstract 

We investigate a forest-fire model with the density of empty sites as con- 
trol parameter. The model exhibits three phases, separated by one first-order 
phase transition and one 'mixed' phase transition which shows critical be- 
havior on only one side and hysteresis. The critical behavior is found to be 
that of the self-organized critical forest-fire model [B. Drossel and F. Schwabl, 
Phys. Rev. Lett. 69, 1629 (1992)], whereas in the adjacent phase one finds 
the spiral waves of the Bak et al. forest-fire model [P. Bak, K. Chen and 
C. Tang, Phys. Lett. A 147, 297 (1990)]. In the third phase one observes 
clustering of trees with the fire burning at the edges of the clusters. The re- 
lation between the density distribution in the spiral state and the percolation 
threshold is explained and the implications for stationary states with spiral 
waves in arbitrary excitable systems are discussed. Furthermore, we comment 
on the possibility of mapping self-organized critical systems onto 'ordinary' 
critical systems. 

PACS numbers: 05.40.+j, 05.70.Jk, 05.70.Ln 



1 



Typeset using REVTgX 



I. INTRODUCTION 



In 1990, Bak et al. introduced a simple model for the spreading of a fire in a forest or 
the spreading of disease in a population ffl. The individuals (sites on a square lattice in 2 
dimensions) can be in one of three states: tree (healthy, excitable), tree on fire (infected, 
excited), and ashes or empty site (immune or dead, refractory). New individuals are 'fed' 
into the system with a small rate p. Whether the third state is considered as death of an 
individual and p consequently as the birth rate of new individuals, or as a state of immunity 
and 1/p as time scale of the loss of that immunity, is a matter of interpretation. In the 
following, we will use the terms tree, fire, empty site and refer to p as the tree growth rate. 
The exact rules of the model [0 were as follows: (i) at each time step trees grow at empty 
sites with a small probability p, (ii) trees on fire will burn down the next time step and turn 
to empty sites, (iii) the fire on a site will spread to the trees at its nearest neighbor sites 
at the next time step. Although originally claimed to be critical in the limit p — > 0, the 
simulations in [0,|J showed that the model does not display criticality. Instead, one could 
observe quasi-deterministic spiral waves of fires. 

In 1992, Drossel and Schwabl introduced the self-organized critical forest-fire model (SOC 
FFM) [|J with the additional rule: (iv) if no nearest neighbor is burning, a tree catches fire 
with a small 'lightning' probability /. Under the condition of a double separation of time 
scales (time between two lightning strokes 1/f 3> time scale of tree growth 1/p > time 
needed to burn down large tree clusters) the model shows critical behavior over a wide 
range of parameter values. The properties of this model were investigated in e.g. [fjUI0,ll|. 

In this paper, we shall investigate a model with the same type of interactions while 
keeping constant the number of empty sites or immune individuals. Their density is the 
control parameter of the model. Parameters like p or / do not enter the model. We will 
show that the model exhibits three phases, two of which can be shown to display the same 
behavior as the abovementioned models. The third phase shows clustering of trees with the 
fire burning at the edges of the clusters. 
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In the 'spiral wave' phase, which exists not only in this model, but in a large number 



of excitable systems (for a review on excitable systems see e.g. PJIH[]), we will point out an 
interesting relation to the non-equilibrium percolation model of JIT|| . There, the tree density 
was the control parameter and the following rules were iterated : (i) Lightning strikes an 
arbitrary site in the system. If the site is occupied, the whole cluster of s trees which 
is connected to this site (by nearest-neighbor coupling) burns down, i.e. the trees of this 
cluster turn to empty sites, (ii) Then s new trees are grown at randomly chosen empty sites 
(including the ones that have just turned empty). The close relation between this model 
and the model treated in the present paper will enable us to identify the mechanism that 
determines the density distribution of excitable constituents and, in particular, the density 
immediately in front of excitation fronts in the spiral-wave phase of excitable systems. 

Furthermore, the model presented in this paper, as well as the model in ||11|| , are examples 
of reformulations of a SOC model in terms of a control parameter (the density of trees or 
empty sites, respectively). Both models indicate that the claim in ||12[ that all SOC models 
can be mapped onto ordinary critical systems exhibiting a sub-critical phase, a critical 
phase with a smoothly varying order parameter and a critical point that separates the two 
phases, is not true in general. Non-equilibrium systems and their phase transitions show a 
much richer behavior than equilibrium systems, with many features that are unknown in 
equilibrium. 



II. THE MODEL 

The model is defined on a d-dimensional hyper-cubic lattice with L d sites. If not stated 
otherwise, we choose d = 2 and periodic boundary conditions in the following. p e L d sites 
are randomly chosen to be empty. The density of empty sites p e is the control parameter of 
the model. The remaining sites are randomly filled with trees and fires. Their densities are 
denoted by p t and pf. It is always p t + p e + Pf = 1- The exact values of p t and pf in the 
initial state do not affect the stationary state, except in the vicinity of certain points which 
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will be discussed later. 

The system is iterated as follows: (i) all trees on fire will burn down the next time step, 

(ii) the fire on a site will spread to the trees at its nearest neighbor sites in the next time step, 

(iii) after each time step the same number of trees that have burnt down grow at randomly 
chosen empty sites (including the ones which have just become empty), thereby keeping p e 
fixed, and (iv) if the fire dies out, a randomly chosen tree catches fire spontaneously 

The motivation of rule (iv) is the following: We want to investigate the system under 
the action of a vanishingly small lightning probability . Since the process which is described 
by rule (iv) may then take very long in real time, but nevertheless can be simulated in one 
iteration step, rule (iv) represents an acceleration of the real process. When calculating 
temporal averages of the fire density, this point has to be considered. 

The reason for choosing the density of empty sites as parameter and not the density of 
trees is the following: Consider a system that consists only of trees and replace one tree by 
a fire. In the next step we have four fires, but only one empty site to grow new trees. This 
extreme example shows that there might be situations in which it is not possible to keep the 
density of trees constant. The density of empty sites, in turn, can always be kept constant 
for arbitrary values within the interval [0,1]. 

In the following sections, we discuss the properties of the stationary state as function of 
the density of empty sites p e . In Sec. |TTT|, we start with a high density of empty sites and 
investigate the region of vanishing fire density. In Sec. [IV|, we lower the density of empty 
sites which leads to a state with spiral waves. A detailed description of the mechanism which 
determines the density distribution in the spiral state is given in Sec. |V|. Sec. [VI| investigates 
the more homogeneous, 'mixed' phase that can be observed after further decreasing the 
density of empty sites. Up to here, only two-dimensional square lattices are considered. 
Sec. |V1I| treats other dimensions and lattice types. In Sec. |V111| , we comment on the issue 
of mapping self-organized criticality onto ordinary criticality. Finally, in the appendix, we 
present some general properties of the order parameter curve. 
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III. REGION OF VANISHING FIRE DENSITY AND CRITICAL POINT 



For p e < 1, there exist only very small tree clusters, and, consequently, if one starts a 
fire by setting on fire a randomly chosen tree, it soon dies out, and one has to start a new 
one. The average number of trees s destroyed by a fire therefore is finite and small, and 
in the thermodynamic limit L —>■ oo the fire density pj equals zero (taking into account 
that rule (iv) is an acceleration of a process which takes infinitely long in real time). In 
our simulations the maximum system size was 6000 2 . With decreasing p e , s increases, but 
still remains finite, and therefore pf = 0. If p e is decreased further, we finally arrive at a 
critical density p^ 1 ~ 59.2% (p^' 1 = 1 — p^ 1 m 40.8%), where s diverges with a power law 
s oc (p e — Pg' 1 ) -5 , with some exponent 5. A snapshot of the system in the vicinity of p^ 1 is 
shown in Fig. [I]. 

The critical behavior close to p^ 1 can be described by exponents which are defined as in 
percolation theory fL3f . The size distribution of tree clusters is n(s) oc s~ T C(s/s max ) with 
a cutoff function C. s max is the size of the largest cluster in the system and diverges for 
p e — > Pg' 1 . The fractal dimension p of the clusters is defined by R(s) oc s 1 '^, where R is 
the radius of gyration of a cluster. The correlation length is given by £ oc (p e — Pe' 1 ) - " '■ 
More exponents can be defined and scaling relations between them can be derived (see e.g. 

The critical exponents found in the simulations (r = 2.14(4), p = 1.95(2), and v = 0.28) 
are the same as in the SOC FFM, when appropriately redefined (for v one has to change 
variables from f /p to p e via p e — p^ 1 oc (f/p) 1 ^ 8 (see |§)). Also the critical density pi' 1 ~ 
40.8% remains the same. This model displays exactly the same critical behavior as the SOC 



FFM and the model of ||11|| , because in a system which is much larger than the correlation 
length, neither the difference between a globally conserved density p e (this model and |TTJ) 
and a density p e which is only conserved on an average (SOC FFM) nor the difference 
between instantaneous re-growth of trees (this model) or delayed re-growth (SOC FFM and 



IIJ) can be seen on length scales comparable to the correlation length. 
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In the stationary state of the SOC FFM s = pp e / f(l — p e ), since in one time step there 
are PtL d f lightning strokes and p e L d p new trees are growing. Therefore, if we measure in 
our model s for a certain value of the control parameter p e we know that its behavior is that 
of the SOC FFM for f/p = p e /s(l - p e ). 



IV. REGION OF FINITE FIRE DENSITY I: SPIRALS 

If we decrease the number of empty sites beyond the critical point p^ 1 , the size of the 
largest forest cluster diverges, and we would expect the fire not to be extinguished any 
more. Therefore, one might expect the fire density pf to behave as an order parameter 
that sets in at Pg' 1 and grows smoothly from zero to finite values, obeying some power law 
Pf oc (p^ 1 — PeY ■ Instead, the following behavior can be observed in the simulations: The 
system restructures, fires gathering to form spirals and the regions of different densities (see 
Fig. [I]) vanishing in favor of a smooth density distribution between the spiral arms (see 
Fig. |2|). The behavior of the spirals is quasi-deterministic. Immediately in front of the fire 
fronts the tree density is very high and immediately behind them it is very low. The density 
distribution is treated in more detail in Sec. [V|. The distance A between two spiral arms is 
finite and constant throughout the system. The fire density for p e = p^' 1 — is also finite. 
Decreasing p e further, the distance between the spiral arms becomes smaller and the fire 
density becomes larger, since it is now easier for the fire to survive. At p e = p^' 3 ~ 54.2%, 
the spiral state breaks down and another restructuring to a new phase takes place, which 



will be treated in Sec. [VI . 

For p e < p^ 1 the fire is able to sustain itself and it is no longer necessary to keep it alive 
by setting on fire randomly chosen trees each time the fire has died out. The 'external field' 
can be set zero. In terms of tree growth and lightning, we have a model with tree growth 
rate p, but without lightning rate /. Therefore, the behavior in the spiral state at some 
density p e is exactly the same as in the forest-fire model of Bak et al. in [|1|] for a certain 
value of the tree growth probability p, since in the thermodynamic limit it does not matter 
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which of the variables is kept fixed. The fluctuations of the densities in the Bak model 
vanish for L — > oo, as do the fluctuations of the number of new grown trees per time step 
in this model. 

In the Bak model, one observes a diverging distance between spiral arms for p — > 0. If 
we start in our model with a spiral state and increase p e beyond p^ 1 (corresponding to a 
decrease of p), the system does not undergo a phase transition at p e = p^ 1 , but chooses to 
retain the spiral state up to p e = p^ 2 ~ 60.8%. This value of p e with A = oo corresponds to 
p = in the Bak model. As function of p the order parameter follows the power law pf oc p 
(explanation in Sec. |V|). 

As already found in , the spiral state does not display criticality or scale invariance in 
the sense of clusters or events on all length scales. For a particular value of p the model does 
not show events on all scales, but only on the scale 1/p, and behaves essentially deterministic. 
However, the model allows similarity transformations, since a system with parameter p\ can 
be obtained by rescaling a system with parameter p 2 . 

What we have found is a discontinuous phase transition which is neither of first nor 
of second order. We find hysteresis (the state of the system for a density in the interval 
[Pe 1 ! Pe' 2 } depends on from which side one approaches this region), an order parameter that 
increases smoothly from zero at p^ 2 to finite values, and critical behavior on only one side 
of the transition (p e > p^ 1 , the side with vanishing order parameter). The order parameter 
curve for a two-dimensional system is shown in Fig. |3|. 

The reason why the critical point p e = p^ 1 is not an ordinary critical point is the same 
TTJ| . From Fig. [I], it can be seen that, in addition to large tree clusters, there exist also 
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large clusters of empty sites. In contrast to ordinary critical phenomena (e.g. percolation), 
there is no homogeneously distributed set of large clusters that could join at p e = p^ 1 to 
form an infinite cluster that spans the whole system. Rather, the largest cluster has to 
compete with all other regions with different densities (which, due to the very nature of 
the dynamics, are nothing more than different growth stages of the largest cluster itself) 
for space in the system. Or, to put it differently, an infinite cluster like in percolation is 



impossible, because the system has to provide space not only for the infinite cluster, but 
also for a large number of 'younger' copies of it. The violation of the hyperscaling relation 
d = h(t— 1) also indicates the inhomogeneous distribution of density in the forest-fire model. 
As explained in ||, the violation is equivalent to the statement that not every part of the 
system contains a spanning cluster at criticality. 



V. DENSITY DISTRIBUTION IN THE SPIRAL STATE 

In order to understand better the density distribution in the spiral state, we shall now 
investigate a state with a single front propagating through a square system with periodic 
boundary conditions (see left side of Fig. |). Since this can be considered to be a section 
of a spiral state (see right side of Fig. f|), and the spiral state, in turn, can be completely 
covered by such sections (if the edge length is chosen to be equal to the distance A between 
two successive arms of the spirals), it is sufficient to understand the density distribution in 
this part of the whole system. Instead of successive fire fronts passing through the small 
section of the system, one might also think of a single fire front which repeatedly leaves the 
section at one end and reenters at the other end due to periodic boundary conditions. 

The tree density immediately in front of the fire front pj? efore has to be larger than or 
equal to the percolation threshold for site percolation on a square lattice p c ~ 0.59, for 
otherwise there wouldn't be propagation. With increasing distance from the fire front, the 
tree density smoothly decreases, until one finally arrives again on the other side of the 
front, where the density takes its lowest value pf ter . If pj? eforo was exactly the percolation 
threshold, the fire would burn down a vanishing fraction of all trees, and the propagation 
speed of the front would be zero. The densities in front of and behind the fire front would 
be equal (pt> efore _ p c _ p^ fter ) ; as claimed in However, the simulations clearly show 
that this is not the case. pj? efore is higher than the percolation threshold, the density pf ter 
in a region after the fire has passed through it is very low, but not equal to zero, and the 
propagation speed of the fire fronts t>fi re is non-vanishing (0.7 ± 0.05 sites per iteration step 



S 



near Pg 2 ). What is the mechanism that determines the density distribution and the values 
of the density immediately in front of and behind the fire front? 

We shall approximate the periodic 'single front' state even further by a coarsened state 
with a finite but large number n of stripes of equal width and different densities parallel 
to the fire front. Each stripe has homogeneous density. Let pj,...,p™ be the densities 
of the n stripes, starting with the highest density. p\ is then the average density of the 
original smooth 'single front' system in the region covered by the ith stripe (p\ pjpefore 
and p™ pf tCT )- We also coarsen time, in that we consider the propagation of the fire front 
from the bottom of one stripe to its top as one time step. Growing of new trees then takes 
place after each (coarsened) time step. If one sets on fire the baseline of the stripe with 
highest density, the infinite cluster in that stripe will surely be set on fire, since it has many 
connections with the baseline, and the propagation of the fire front effectively causes the 
removal of the infinite cluster from this stripe. 

The propagation of the fire front therefore can be modeled as follows. First, we identify 
the infinite cluster in the stripe with highest density and remove it from the system. Doing 
this one has to respect the boundaries of the stripes, although the infinite cluster of course 
extends into the neighboring stripes. The strength of the infinite cluster at density p\ is 
denoted P(p])- (We now define also the density p™ +1 of the stripe which contained the 
infinite cluster after its removal, i.e. p" +1 = p\ —P(p\).) Second, the P(pj)A 2 /n trees of the 
infinite cluster are redistributed randomly amongst the empty sites of the whole system of 
size A 2 . If the system is to be stationary, the stripes thereby just exchange their densities, 
i.e. the stripe with second highest density p\ now assumes the highest density p], and so 
on. From that condition, one can easily derive the equation 



for i — 2, . . . , n + 1. The last factor on the r.h.s. represents the fraction of trees of the 
infinite cluster that are re-grown in the stripe with density p\. We finally obtain 



Pi 1 = Pl + (Pi - P n t +1 ) ■ (1 - P\)/{n{l ~ Pt) + P \ ~ P n t +1 ) 



l-p\ _1-P 2 t 




(1) 



i-P 2 i- Pi 
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Together with p t = (1/n) Y%=i Pt we have n equations for n + 1 densities. These equations 
were already derived in ]TI| for a related model. The average density in a system with n 
stripes is then: 



Pt 



Sir 



1 p\ - Pt +1 



n(((l-pi)/(l- P r +1 ))- 1/n -l) 
In the system that we are really interested in, the density varies smoothly, so we have to 
consider large n. In it was argued that p\ has always to be greater than or equal to a 



certain constant pi « 0.625 > p c . If this were not the case, a traversing fire front would leave 
behind large tree clusters which would lead to inhomogeneities that prevented the next fire 
front from passing through in the same quasi-deterministic way than the first one. Likewise, 
has always to be greater than or equal to another constant pf = p* t — P(p* t ) ~ 0.078. 
p\ is always larger than the corresponding percolation threshold p c . These constants have 
a more fundamental significance independently of this model and of any states with spiral 
waves (for details see 0,0). The case p) = p* and p™ +1 = pf 3 corresponds to the lowest 
possible overall density and therefore to the state with infinite spirals at p = or p e = p ^ 2 . 
In this case, we have p^ efore = p* and pf tcr = pf 2 . 

The overall density in the spiral state can then be written as 

^before „aftcr 

_ -i Pt Pt 

Pt J- i „bcforc , 



hm n ^ 00 n((^^)- 1 /n_i) 



before after 

Pt Pt (2) 



lim^oon^l - ln( _' ltcI ) n ± . . . - 1 

1 ft 



^before ^.after 

I — L ™ 



after 
t ' 
before 



Mr - 

For p^ eioie = p* and pf ter = p^°, this density should be equal to the density pt' 2 of the state 



with spirals of infinite extension. With the values of p% and pf measured in [JTTJ] (0.625 
and 0.078 for the square lattice, 0.533 and 0.062 for the triangular lattice), one arrives at 
^.calculated = q 392 ± q.002 for the square lattice and pf - calculated = 0.325 ± 0.002 for the 
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triangular lattice, in excellent agreement with the values p£' 2 = 0.392 ± 0.002 for the square 
lattice and pi' 2 = 0.323 ±0.005 for the triangular lattice, measured in the spiral state of this 
model. 

Since Eq. (|I|) represents a geometric series, the density as function of the distance from 
the fire front is given by 

(i before \ e/A 

with Pt (0) = pf tCT > p t °° and p t (A) = p^ cfore > p*. 

With our new knowledge about the nature of the spiral state, we can derive an equation 
for the distance A between the fire fronts as function of the tree growth probability p. Let 
the speed of the fire fronts v^ re (p^ eiore ) be measured in sites per iteration step. Like P it 
depends on the density in front of the fire front p^ efore . The amount of matter burnt in unit 
time is then P(pj? efore )vfi re (pj? eforc )A. This has to be equal to the number of growing trees 
per unit time (1 — p t )pA 2 , leading to the relationship 

= P(pr° re Kre(pr° rC ) 1 m 

I-Pt v u 

This confirms the observed scaling behavior A oc p -1 and additionally delivers the constant 
of proportionality, being 0.63 ±0.05 for p^ efore = p* and in good agreement with simulations. 
Since the fire density is inversely proportional to the distance between the spiral arms, we 
find pj oc p, i.e. the order parameter exponent j3 equals one. This can also be seen from the 



equality p/ = pp e jlql , which states that the number of burnt trees has to be equal to the 
number of new grown trees in the stationary state. 

If one chooses a fixed site behind the front and wants to know the tree density in this 
region as function of the time t until at T = A/vr tc the next front passes through, one has 
to replace in Eq. (||) x by v^ re t and arrives with Eq. (|]) at 

r / u r V 1 P* 

/ -| ^before \ pf 

Pt (t) = i-(i- P fn 



i „after ; 
1 ~ Pt ) 

which with Eq. (|) and P(p t bcforc ) = p|? cfore - pf ter leads to 
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p t (t) = 1 - (1 - pf e *)e- pt . (5) 

Eq. (|) can also be derived easily from dp t /dt = p{l — p t ), first stated in 0. This shows 
that our picture is in accordance with the basic equations describing the spiral state. 

One additional point concerning Eq. (Q) should be mentioned. Eq. (Q) relates A, p t , and p. 
If one regards e.g. A as fixed, only the product (1 — p t )p is determined. This reflects the fact 
that in a 'single front' state with fixed edge length A one can have different overall densities 
Pt- The other quantities adjust themselves automatically. However, in the 'real' spiral state, 
fixing A determines invariably all other quantities. Since the spiral state can be constructed 
by connecting 'single front' systems, this seems to be a contradiction. The solution of this 
apparent contradiction is that one has neglected the spiral centers. They yield a second, 
unfortunately unknown, relation between A, p t , and p. This relation is brought about by 
the rotation of the spiral centers, pt and p determine the angular velocity u of the rotation, 
which, in turn, determines A via A = v^ rc T = v^ vc 2ti / uj . For the calculation of pi' 2 with 
Eq. (H), using the decomposition into 'single front' states, it is justified to neglect the spiral 
centers, because in the thermodynamic limit their number density is zero. Nevertheless, 
one has to be aware of the fact that, although not important for calculating p^' 2 , the spiral 
centers are the 'pacemakers' of the spirals and therefore responsible for the magnitude of A. 

Another interesting point is that the density pf 2 can also be found from an extremum 
principle first stated in Q. There, the extremum principle was erroneously used to determine 
the critical density p c t ' 1 of the SOC FFM. The fire was believed to destroy as much trees as 
it can at the critical point, but the result ~ 39% was in contradiction to the measured value 
p c t ' 1 f=a 40.8% |5]|||^],|3|]. However, the equations derived in [[| from the extremum principle 
can easily be seen to be equivalent to Eq. (§) for the spiral state. Therefore, the principle 
yields the correct critical density, but for a different, then unknown, phase of the model. 

From the results found in this section one can draw some conclusions for excitable media 
in general. In many excitable systems stationary states with spiral waves can be found. 



Famous examples are e.g. the Belousov-Zabhotinski (BZ) reaction or the electro 
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physiological activity of heart tissue ||17|| . If the spirals are to sustain themselves in a 
stationary state, the density of the excitable constituents (the 'fuel') in a region imme- 
diately before a fire front passes through it pj? efore has always to be larger than or equal 
to some threshold p^T, which, in turn, is larger than the percolation threshold p c for that 
particular situation. The percolation threshold can in principle be measured by preparing 
a homogeneous system with a certain density of excitable constituents and 'exciting' one 
edge of the system. If the resulting excitation front dies out before reaching the other end 
of the system, we are still below the percolation threshold, and if it reaches the other end in 
finite time we are above. At the percolation threshold the front barely survives and needs an 
infinitely long time to reach the other end. The overall density of the excitable constituents 
has to be above the value p£' 2 , below which no excitation can be sustained. 

A system that displays spiral waves or avalanches (concentric waves, target patterns) 
depending on the density of excitable constituents are e.g. populations of Dictyostelium 
Discoideum amoebae ||18||. There, circular waves of signaling activity emanating from pace- 



makers are found for low cell densities, whereas for high cell densities one observes spiral 



waves. 



VI. REGION OF FINITE FIRE DENSITY II: THE MIXED PHASE 

With decreasing p e , the distance between the spiral arms becomes smaller and the spirals 
finally break up to single fronts. This change occurs continuously. The fronts, which are 
more irregular than the spirals, have also been observed in the model of Bak et al. De- 
creasing p e further, the coherence length of the fronts soon becomes comparable with the 
lattice constant, and the system reorganizes itself discontinuously into another state, which 
constitutes the third phase of this model. This happens at p e = p^ 3 ~ 54.2%. The simula- 
tions show that the trees in the new state tend to 'cluster' (to form regions with higher tree 
density) with the fire burning at their edges (see Fig. ^). To sustain this type of structure, 
a certain minimum density of fires is needed. At p e = p^ 3 , pj jumps from 2% to 10%. If 
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one lowers p e further, the size of the 'clusters' decreases, and the fire density increases. The 
system as a whole becomes more homogeneous (trees, fires and empty sites are 'mixed'). For 
p e = 0, we observe pt — pf — 1/2 for all dimensions and lattice types. The order parameter 
curve starts linearly at p e = with a slope that depends only on the number of nearest 
neighbors z (explanations in the appendix). If we start with small p e and traverse the phase 
transition in the other direction, it takes place at different p e (p^ 4 ~ 54.7% > pi' 3 ), i.e. one 
has a first-order phase transition with hysteresis. The order parameter as function of the 
density of empty sites is shown in Fig. §. 

VII. DIMENSIONS OTHER THAN TWO AND DIFFERENT LATTICE TYPES 

Since in one dimension there cannot exist spiral waves, we expect a simpler phase diagram 
than in two dimensions. For p e = we are dealing with a completely deterministic one- 
dimensional cellular automaton where each site can be in one of two states (automaton no. 
54, according to the classification scheme of Wolfram [|19|]). In the stationary state each tree 
has at least one fire as neighbor and vice versa. Strings of more than two trees or fires are not 
stable, since sooner or later they would be invaded by fires. The stationary state is periodic 
with period two. The introduction of empty sites makes the automaton non-deterministic, 
because the trees can now choose where to grow. This variant has not been investigated 
so far. In the simulations, system sizes of up to L = 10 7 were used. The measured order 
parameter curves for the one-dimensional case can be seen in Fig. [| In addition to the usual 
nearest-neighbor interaction, we simulated also a variant where the fire is allowed to jump 
over one empty site if necessary. The curves start linearly at pf = pt = 1/2 with different 
slopes (explanation see Appendix). For increasing p e one observes the same phenomena as 
in two dimensions. The tree clusters (strings in ID) become larger and are accompanied by 
fire on at least one of their ends. At p e = p c e ' A ~ 20% (pl' A ~ 26% for next-nearest neighbor 
interaction) the density of empty sites is too high for this structure to survive, and the fire 
density drops from a finite value of approximately 20% to zero, since the system does not 
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have the possibility to rearrange itself into a spiral wave phase. Instead, now isolated chunks 
of trees are ignited and burn down in finite time. In the reverse direction, pf remains zero 
until p e = p^ 1 = 0, i.e. the hysteresis in one dimension is maximum. In the vicinity of 



the point (p e = 0, pf = 0), the critical behavior of the SOC FFM in one dimension |20| is 
reproduced. Since there exists no spiral phase in one dimension, there are no pf 2 and pf 3 . 

The two-dimensional simulations were also done with a triangular lattice and yielded the 
same behavior as for the square lattice (for the order parameter curve see Fig. ||). The values 
of pg' 1 ' 2,3 ' 4 and the slope of the order parameter curve for p e — > of course are different (see 
Tab. |). 

In three dimensions, we could identify a sub-critical phase as well as a mixed phase. 
The order parameter curve is plotted in Fig. [7|. Due to small L (L < 300) the critical 
exponents in the subcritical phase and the critical density pf 1 could not be measured with 
sufficient accuracy. The mixed phase in three dimensions shows analogous behavior to the 
two-dimensional mixed phase and does not display new phenomena. The fire density curve 
could not be measured for very small pj (pj < 0.2%), due to finite size effects. If one 
extrapolates the fire density curve to the point pf = 0, one arrives at 77.3%, therefore, pf 4, 
lies between the last simulated density 76.9% and 77.3%, i.e. pf 4, = 77.1% ±0.2%. The fact 
that there exists a gap between the end point of the critical phase pf 1 ~ 78.1% (taken from 
||) and the end point of the mixed phase p^' 4 leaves open the possibility of the existence of 
a third phase containing the three-dimensional analogon of spiral waves (scroll waves). Due 
to the finite size of our sample, however, we could only observe two-dimensional, flat fire 
fronts in that region. 

The results for all simulated lattices and dimensions can be seen in Tab. |. 

VIII. MAPPING SELF-ORGANIZED CRITIC ALITY ONTO CRITIC ALITY 

The model we have investigated in this paper is an example of a system which is far from 
equilibrium. It shows a wealth of interesting structures depending on the density of empty 
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sites p e . For large p e , we find a region of vanishing fire density which contains the critical 
behavior of the SOC forest-fire model. For lower density of empty sites we obtain the spiral 
waves of the Bak et al. model. Finally, for even lower p e , we find a phase in which the 
trees show the tendency to form clusters with the fire burning at their edges. The transition 
between the second and the third phase is of first order, whereas the transition between 
the first and the second phase is rather unconventional. We find hysteresis and critical 
behavior on only one side of the transition. The close relation of the spiral wave phase with 
the synchronized phase of the model of JO] allows us to understand the mechanism which 
determines the density distribution in the spiral state of not only this model, but of arbitrary 
excitable systems. In particular, it yields the density in front of the excitation fronts, the 
minimum density of excitable constituents that is necessary to sustain the excitation, and 
the factor of proportionality between the distance of the spiral arms and the 'tree growth' 
rate p. 

Apart from being interesting in its own right, the model presented in this paper (like 
the model of [IT]) has to be seen also in the context of the general claim of [12], that the 



critical points of SOC models can be regarded as ordinary critical points of second order 
phase transitions. It was claimed in [T2J that this should be possible for all SOC models. 



There were also given instructions on how to achieve this goal for some particular models 
including the forest-fire model. 



The results in this paper and in ||11|| , however, disprove this hypothesis. While the 
sub-critical side of the transition in both models behaves as expected, the other side (the 
side with a supposed non- vanishing order parameter) always exhibits surprising features. In 



the model of [[11]] we found a whole critical region and lots of first-order phase transitions, 



whereas the actual model shows hysteresis and no criticality at all. 

Although these 'negative' results can not strictly rule out the possibility that for some 
further slight change of the rules one might succeed in obtaining the usual 'decent' behavior 
of the order parameter, the argument at the end of Sec. [TV| , together with the simulation 
results of this paper and of |ll[ make it seem very unlikely. We consider it to be more 
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probable that a subset of SOC models are 'only' ordinary critical models 'in disguise'. With 
these models the mapping proposed in |]12[ should be possible (the results in [2TJ also seem to 
point in that direction), while there also exist genuine SOC models for which this procedure 
can not be carried through. The forest-fire model seems to belong to the second class. 

We suggest that similar phenomena and difficulties (critical regions, hysteresis at the 
critical point, critical behavior on only one side) will also be found in other models of SOC, 
as e.g. the sandpile, earthquake, and evolution models. 

This work was supported by the Deutsche Forschungsgemeinschaft (DFG) under Con- 
tract No Schw 348/7-1. 

We thank B. Drossel for fruitful discussions. 



IX. APPENDIX: SOME PROPERTIES OF THE ORDER PARAMETER CURVE 

All order parameter curves considered in the previous sections started at the point P/ = 
p t — 1/2, p e — 0. This starting point is independent of dimension or lattice type, which can 
be seen as follows. For p e = each tree will sooner or later be set on fire by one of its z 
nearest neighbors. Once having set on fire the tree, the fire jumps forever between this site 
and its neighbors, because its neighbors, after burning down, have to become trees again 
in the next time step (there are no empty sites) . The dynamics of these 1 + z sites then is 
fixed. After each time step trees and fires change places. But the fire can still propagate to 
other sites, so that finally the whole system consists of such 'blinking' regions. The state is 
periodic with period 2 and the average fire and tree densities have to be 1/2. For a random 
initial state the densities are 1/2 also without averaging. 

The slope of p/(p e ) for p e — > depends only on the number of nearest neighbors z and lies 
always within the interval [—1; —1/2]. This can be seen as follows. If one takes a stationary 
state at p e = and inserts some empty sites by removing an equal number of trees and fires, 
the system will not remain in this state (with a slope of —1/2), but will adjust itself to a 
new stationary state with an even lower number of fires, since the spreading conditions for 
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the fire are now worse than before due to the inserted empty sites. The magnitude of this 
effect depends only on the coordination number z. The larger z, the smaller is the effect. 
On an infinite-dimensional lattice, the fire density will not readjust itself at all, because no 
fire can feel the empty sites. The slope in this case is the maximum possible slope —1/2. 
This can also be seen mathematically. All fires have been trees in the preceding time step, 
therefore pf < p t in the stationary state. With pf + p e + p t = 1 it follows pf < 1/2 — p e /2. 
Since for p e = both sides of the last inequality are identical, one can differentiate and finds 
dpf / dp e \ Pe= o < — 1/2. 1/2 is the upper bound for the slope. 

The lower bound for the slope equals —1. This value is assumed for one dimension, 
as will be shown in the following by considering the effect of introducing one empty site 
into the stationary ID state at p e = 0. Since the one- dimensional lattice has the smallest 



possible number of neighbors, its slope represents the lower bound. As argued in Sec. [VI H , 
in the one-dimensional stationary state at p e = no more than two neighboring sites can be 
in the same state. In a symmetric state, where both neighbors of a tree are fires and vice 
versa, i.e. where fires and trees are strictly alternating, an empty site does not modify its 
neighborhood. Since the effects of the empty site are only local, it is sufficient to consider 
only three cases: One pair of sites in the same state (see left part of Fig. |S[), collision of two 
pairs (see middle part of Fig. §), an d annihilation of two pairs (see right part of Fig. ||). 
The possibilities of placing the empty site in these cases can be further reduced to the cases 
shown in the second line of Fig. || The reason is that, since in each time step n trees burn 
down and are refilled into n + 1 empty sites, with probability ~ 1 for large n the empty 
site at time t was a fire at time t — 1, i.e. at time t the empty site occupies a place which 
would have been a tree in the unperturbed state. Therefore, we place the empty site at a 
tree site in the second line of Fig. |8|. If instead we had chosen the neighboring tree sites, we 
would only have generated mirror-symmetric forms of the configurations shown. In the third 
line the empty site has disappeared from our small section of the system and has become 
a tree with probability ~ 1 for n»l. In the left case, the pair of sites moves one lattice 
spacing due to the presence of the empty site. In the middle case, two pairs which are 



separated by only one lattice site collide due to an empty site. In the right case, an empty 
site causes two neighboring pairs (which may have collided earlier) to annihilate each other 
and thereby symmetrize the region around the empty site. Therefore, sooner or later all 
pairs will have annihilated, and the entire system (with the exception of the single empty 
site) is in the symmetrical state. This evolution to the symmetrized state was also observed 
in the simulations. 

We can therefore restrict ourselves to consider the consequences of the introduction of 
empty sites into such a highly symmetrical state. This has the advantage that there is 
no need to consider the neighboring sites, because they are not affected. There are two 
possibilities to introduce an empty site. The first one is to remove a tree. In the next 
time step this tree would have become a fire which it cannot do now. Instead, with very 
high probability a new tree is grown at this empty site. During the re-growth of trees the 
empty site effectively moves to another site in the system (which would have been a 'tree 
site'). While there is one excess tree in one part of the system, there is one tree less in 
another part of the system, so the total number of trees does not change. The number of 
fires, however, has decreased by one, since we prevented the 'birth' of one fire. The second 
possibility is to remove a fire initially. If that happens, the empty site will have moved after 
the tree-growth phase to a 'tree site', and we have the first case again. Therefore, in each 
case the introduction of empty sites takes place completely at the expense of fire sites and 
dp f /dp e \ Pe=0 = -1. 

The measured results for the slope of the order parameter curve are shown in Tab. |. One 
can see that the slopes of the two-dimensional triangular lattice and the three-dimensional 
hyper-cubic lattice are the same, since they have the same coordination number z = 6. 

As mentioned earlier, not all tree densities between and 1 are possible in the stationary 
state. If one transforms the p/(p e ) curves to Pf(pt) curves via pj + p t + p e = 1, one can read 
off the diagram the allowed p t values (see Fig. |9|). For these values it is possible to simulate 
the model with fixed p t . The results do not differ from the results with fixed p e . Also in 
this case, one has to be careful with the initial conditions, since there are sometimes two 
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stationary states for the same value of p t , and the initial state determines which of the two 
possible stationary states will be chosen by the system. 
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TABLES 



lattice type ID ID (*) 2D 2D triangular 3D 

p c / 0% 0% 59.2(1)% 66.4(4)%! 78.1(1) 

p c / - - 60.8(5)% 67.7(5)% 



c / - - 54.2% 61.7% 



p c / 20% 26% 54.7% 62.2% 77.1(2)% 

1.00(1) 0.85(1) 0.67(1) 0.54(1) 0.55(1) 



dpj_ 

dpe 



p e =0 



TABLE I. The simulation results for various dimensions and lattice types i 1 = taken from [8]). 
In the case ID (*) the fire was allowed to jump over one empty site, if necessary. 
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FIGURES 



FIG. 1. Snapshot of the stationary state in the 'SOC phase near the critical density 
Pe' 1 ps 59.2%. L = 2000 and p e = 59.8%. Trees are black and empty sites are white. 



FIG. 2. Snapshot of the stationary state in the 'spiral wave' phase for L = 1000 and p e = 59%. 
Trees are grey and empty sites are white. The fires are black, but difficult to see. They are located 
at those lines where the density of trees changes abruptly. 
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FIG. 3. The order parameter fire density pf as function of the density of empty sites p e in two 
dimensions ( — , square lattice; • • •, triangular lattice). The magnified sections in the insets show 
the two transitions of the two-dimensional system on a square lattice. To the left, one can see the 
first order phase transition at p^ 3 and p^ 4 between the mixed phase and the spiral phase, and to 
the right the 'mixed' phase transition at p^ 1 and p^ 2 between the spiral phase and the SOC phase. 
The arrows indicate the directions in which the transitions are traversed. 
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FIG. 4. A spiral state at density p e = 59% with distance A between successive fire fronts can 
be covered by quadratic 'single front' sections of edge length A. Trees are grey, fires are black, and 
empty sites are white. 
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FIG. 5. Snapshot of the stationary state in the 'mixed' phase near the critical density 
Pg' 4 ps 54.7%. L = 180 and p e = 54.6%. Trees are grey, fires are black, and empty sites are 
white. 
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FIG. 6. The order parameter fire density pf as function of the density of empty sites p e in one 



dimension ( — , nearest neighbor interaction; • • •, the fire is allowed to jump over one empty site if 
necessary). For p e > p^ 4 « 0.2 and p e > p^ 4 s=s 0.26, respectively, the fire dies out. In the reverse 



direction, pf remains zero until p e = p°' = 0. 
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FIG. 7. The order parameter fire density pf as function of the density of empty sites p e in three 
dimensions. The extrapolation of the fire curve yields the critical density p 1 ^ 4 ss 77.1% ± 0.2%. 
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FIG. 8. The possible effects of the introduction of one empty site into a one-dimensional system 
at p e = 0. 
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FIG. 9. The order parameter fire density pj for all simulated dimensions and lattice types as 
function of the tree density p t . In the case ID (*) the fire was allowed to jump over one empty 
site, if necessary. 
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